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To  determine  whether  a  flawed  solid  propellant  would  be  safe  for  operational  use,  one 
must  first  verify  whether  the  flaw,  which  can  take  form  as  a  crack  or  fracture  (a  volume 
between  two  propellant  surfaces)  or  a  debond  (a  volume  between  propellant  and  liner),  will 
propagate  or  simply  burn  out  after  it  ignites.  A  prediction  of  whether  the  flaw  will 
propagate,  and  to  what  extent  it  will  propagate,  can  be  made  by  calculating  the  pressure 
distribution  inside  the  burning  flaw  and  then  the  resultant  stress/strain  field  generated  in  the 
solid  propellant.  However,  to  precisely  predict  the  pressure  distribution  and  stress/strain 
field,  the  flow  field,  especially  at  the  crack  tip,  has  to  be  thoroughly  understood.  The  work 
described  herein  mainly  focused  on  studying  the  gas  dynamic  behavior  inside  a  simulated 
solid  propellant  flaw  using  a  computational  fluid  dynamics  approach.  In  this  effort,  a  finite - 
volume,  density-based  Navier-Stokes  solver  called  Loci-CHEM  was  used.  The  code 
replicated  experimental  results  with  reasonable  accuracy  and  showed  little  sensitivity  to  grid 
resolution  and  gas  properties  assumptions.  The  pressure  distribution  was  affected  by  the 
assumed  level  of  turbulence  intensity,  which  has  led  to  the  development  of  a  linearized  flow 
instability  model  that  may  be  used  in  the  future  to  better  predict  turbulence  inside  burning 
flaws.  Future  efforts  will  use  fully  coupled  fluid  dynamics  and  structural  dynamics  codes  to 
more  accurately  predict  flaw  propagation. 
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bum  rate  coefficient 
constant  pressure  specific  heat 
turbulent  kinetic  energy 
burn  rate  exponent 
pressure 
grid  spacing 
turbulent  dissipation 
ratio  of  specific  heat 
specific  dissipation 


Nomenclature 


I.  Introduction 

I  hrough  time,  a  solid  propellant  goes  through  slowly-evolving  chemical  changes  (e.g.,  diffusion  of  active 
chemical  species,  energetic  material  decomposition,  etc.)  that  result  in  a  reduction  of  mechanical  properties.  As 
they  age,  a  substantial  portion  of  the  propellant  grain  can  become  either  brittle  or  soft,  depending  upon  its  proximity 
to  liners  and  insulators  and  upon  the  propellant’s  formulation.  Additionally,  cumulative  damage  can  occur  as  the 


1  Aerospace  Engineer,  Propulsion  Engineering  Division,  Member  AIAA. 

2  Director,  Propulsion  Engineering  Division,  Associate  Fellow  AIAA. 

1 

American  Institute  of  Aeronautics  and  Astronautics 
Distribution  A:  Approved  for  public  release;  distribution  unlimited. 


motor  is  subjected  to  rough  handling  loads,  temperature  cycling  and  slumping  from  gravitational  accelerations.  As  a 
result,  it  is  critically  important  to  understand  the  actual  state  of  aged  propellant  prior  to  its  use  in  the  field  or  on  a  test 
stand.  The  use  of  X-rays,  computed  tomography  (CT)  scans,  and  other  non-destructive  evaluation  (NDE)  methods 
is  a  routine  part  of  many  motor  programs.  Those  evaluations  do  not  measure  propellant  mechanical  properties,  but 
they  can  detect  flaws  (e.g.,  propellant  cracks  along  the  bore  and  in  slots,  bondline  de-bonds,  etc.)  that  pose  a  serious 
risk  to  a  motor’s  integrity.  To  do  that,  propellant  samples  are  often  cut  from  aged  motors  so  that  mechanical 
properties  can  be  measured  directly.  These  tandem  approaches  have,  for  the  most  part,  protected  against  failures  of 
aged  motors  and  allowed  questionable  assets  to  be  removed  from  the  inventory. 

Whenever  a  flaw  is  detected  in  a  motor,  the  fundamental  question  that  decision  makers  must  contend  with  is 
whether  the  detected  flaw  will  propagate.  The  logical  follow-on  question  to  that  one,  assuming  the  flaw  will 
propagate,  is  whether  it  will  propagate  to  the  motor’s  destruction.  The  answers  to  those  questions  are  complex, 
incompletely  understood,  and  usually  little  more  than  educated  guesses. 

The  reason  that  is  so  is  because  predicting  flaw  propagation  in  a  propellant  grain  requires  detailed  understanding 
of  at  least  three  distinct  engineering  domains.  For  accurate  prediction,  one  must  understand  the  transient  fluid 
dynamics  and  flame  spreading  behavior  within  the  flaw  volume;  the  resultant  stress,  strain  and  strain  rate  field  that 
develops  in  the  nonlinear  viscoelastic  propellant;  and  the  fracture  mechanics  behavior  at  or  near  the  crack  tip  region, 
which  acts  as  a  stress  concentration  riser.  More  importantly,  for  the  highest  accuracy  prediction,  each  of  these 
domains  must  be  dealt  with  in  a  fully  coupled  manner. 

The  disciplines  have  reached  a  maturity  level  that  can,  with  reasonable  accuracy,  numerically  represent  the 
majority  of  the  physical  phenomenon  associated  with  it.  Thus,  the  overall  goal  driving  the  efforts  summarized  here 
is  to  develop  a  fully-coupled  multi -physics  computational  tool  that  used  the  best  representative  tool  from  each  of  the 
three  domains.  Such  a  tool  would  give  decision  makers  more  information  on  which  to  base  their  go/no-go  decisions 
and  would  do  so  with  a  higher  confidence  factor.  The  first  step  in  creating  that  tool  was  to  verify  that  a  modern 
CFD  tool  could  accurately  simulate  the  combustion  and  flow  inside  a  burning  flaw.  Thus,  the  objective  of  this  study 
was  simply  to  demonstrate  the  feasibility  of  using  the  density-based,  finite-volume,  parallel  Navier-Stokes  solver 
called  Loci-CHEM  code  to  accurately  simulate  the  fluid  dynamic  behavior  inside  a  simulated  propellant-liner 
debond.  A  comparison  was  made  between  the  Loci-CHEM  results  and  a  set  of  experiments  that  was  conducted  in 
the  early  1970’s.  This  validation  was  crucial  for  determining  the  numerical  requirements  (grid  resolution, 
sensitivities  to  assumptions,  etc.)  to  achieve  high  accuracy. 


II.  Description  of  the  Experimental  Data  in  the  Literature 


A  comprehensive  set  of  labscale  experiments  was  carried 
out  in  the  early  1970’s  in  which  the  researchers  developed 
an  experimental  apparatus  to  simulate  the  combustion  gas 
flow  inside  a  propellant-liner  debond  and  measured  the 
pressure  distribution  along  the  simulated  flaw1. 


The  test  article,  reproduced  from  Ref.  1  in  Figure  1,  was 
an  adjustable  wedge  that  had  a  propellant  slab  on  its  lower 
half  and  three  pressure  transducers  installed  on  its  upper 
half.  Two  0.25”  thick  plexiglas  side  pieces  were  used  to 
enable  photographic  visualization  of  burning  processes.  The 
propellant  slab  was  0.1”  thick  (measured  prior  to  ignition) 
and  the  dimensions  of  its  burning  surface  were  3”  wide  and 
4.5”  long.  The  pressure  transducers  were  placed  at  0.19, 

1.72  and  3.72”  from  the  debond  tip.  In  order  to  simulate  the 

externally  imposed  chamber  pressure  that  exists  in  a  rocket  motor,  the  entire  apparatus  was  placed  inside  a  chamber 
that  included  pressure  control  devices  such  that  the  flow  Mach  number  inside  the  simulated  flaw  could  be  controlled 
by  simply  varying  the  wedge  angle.  Both  nitrogen  and  ah'  with  a  range  of  simulated  pressure  were  used  to  pressurize 
the  chamber.  The  debond  angles  and  corresponding  chamber  pressures  used  by  the  researchers  are  summarized  in 
Table  1. 


Figure  1.  Debond  test  article  used  in  Ref.  1 
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To  complement  that  experimental  study,  Jacobs,  et.  al.  also  conducted  a  one -dimensional  steady-state  analysis  of 
a  burning  flaw  using  a  closed-form  analytical  method.  A  comparison  was  made  between  the  numerical  solutions 
and  the  experimental  results.  The  detailed  explanations  of  the  analytical  method  used  can  be  found  in  Ref.  1. 

One  discrepancy  concerns  the  propellant  used  in  [1],  The  authors  reported  using  an  aluminized  double  base 
(DB)  propellant  that  was  described  as  (then)  Thiokol’s  TP-H1011  formulation.  However,  the  Chemical  Propulsion 

Information  Agency’s  (CPI A)  solid  propellant  database 
entry"  for  TP-H-1011  shows  that  it  is  an  aluminized  AP 
(ammonium  perchlorate)  composite  propellant.  Beyond  the 
differences  in  propellant  composition  and  resultant 
combustion  gas  properties,  there  was  a  substantial  difference 
in  burn  rates  reported  in  Ref.  1  and  Ref.  2,  which  obviously 
result  in  differences  between  the  CFD  predicted  pressure 
distribution,  the  experimental  data,  and  the  closed-form 
analytical  solution  shown  in  Ref.  1.  Ultimately,  the  gas 
properties  and  burn  rate  variables  reported  in  Ref.  1  were 
used  in  the  present  study  to  allow  direct  comparison  between 
CFD  results,  the  closed-form  analytical  solutions,  and  the 
Table  1.  Experimental  set  points.  experimental  data 


Wedge  Angle 

Pressures  fp  si) 

Flow  Regime 

4.75° 

170,  300,  400 

Subsonic 

3.17° 

332 

Subsonic 

1.9° 

550,  600 

Subsonic 

0° 

245 

Sonic 

III.  Computational  Tool:  Loci-CHEM 

CHEM  is  a  density-based,  finite-volume,  parallel  Navier-Stokes  solver  for  generalized  grids  using  the  Loci 
framework  developed  at  Mississippi  State  University  in  part  through  NASA  and  NSF  funded  efforts3.  Because  it  can 
use  either  structured  or  unstructured  grids,  Loci-CHEM  is  convenient  for  complex  geometries.  It  is  second-order 
accurate  in  both  space  and  time,  and  uses  high  resolution  approximate  Riemann  solvers  to  solve  turbulent  flows. 
CHEM  can  solve  chemically  reacting  flows,  from  simple  equilibrium  to  complex  multi -species,  multiple-reaction 
finite  rate  chemistry  models.  Several  kinetics  models  are  already  built-in,  and  have  been  used  to  model  the  internal 
flow  of  LOX-hydrogen  propellant  rockets4.  However,  its  capability  for  modeling  the  combustion  and  internal  flow 
of  solid  rocket  motors  is  still  under  development.  There  are  three  main  features  of  modeling  solid  propellant 
combustion  absent  from  the  basic  distribution  of  Loci-CHEM. 

First,  a  boundary  condition  that  models  the  mass  addition  from  the  burning  surface  is  not  available  in  the  basic 
distribution  of  Loci-CHEM.  However,  an  extension  module  was  developed  by  Tetra  Research  Corporation  that  adds 
a  pressure  dependent  mass  flux  boundary  condition  based  on  the  Apn  theory5.  This  module  was  used  in  this  study. 

Second,  the  program  does  not  support  moving  boundaries  and/or  meshes  that  are  needed  to  account  for  the 
propellant  surface  regression.  Another  module  that  will  add  the  capability  is  currently  being  developed,  but  was  not 
ready  for  use  in  this  effort.  Instead,  a  quasi-steady  state  assumption  was  used  where  the  propellant  surface  was 
assumed  to  be  ignited  everywhere  instantaneously.  Commonly  used  in  standard  internal  ballistics  calculations,  the 
quasi-steady  assumption  was  reasonable  given  that  the  summary  data  in  Ref.  1  were  taken  shortly  after  ignition. 
Detailed  photographic  studies  assured  those  researchers  that  the  entire  propellant  surface  was  ignited  quickly.  Also, 
a  special  high  burn  rate  coating  was  applied  to  each  propellant  sample  to  insure  rapid  ignition. 

Lastly,  Loci-CHEM  does  not  support  multi -phase  flow,  when  in  fact  some  of  the  solid  particles  released  from  the 
surface  do  not  burn  instantaneously.  Thus,  in  the  strictest  sense,  the  flow  cannot  be  treated  as  purely  Eulerian. 
However,  it  is  believed  that  non-Eulerian  effects  to  the  overall  pressure  distribution  inside  the  debond  simulator  are 
negligible. 


IV.  Results 


A.  Model  Implementation 

Menter’s  Shear  Stress  Transport  (SST)  model  built  on  a  baseline  model,  namely  the  baseline  (BSL)  model6,  was 
used  in  the  current  effort.  The  model  combines  the  Wilcox  k-co  model  in  the  near-wall  region  with  the  k-s  model  in 
the  outer  part  of  the  boundary  layer.  An  assumption  of  equilibrium  chemistry  and  constant  viscosity  were  used  to 
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model  the  combustion  gas.  The  three-dimensional  mesh  was  generated  using  GRIDGEN  Version  15.  An  anisotropic 
mesh  was  generated  by  extrusion  to  accurately  represent  the  boundary  layer  along  solid  surfaces.  The  remainder  of 
the  computational  domain  was  meshed  using  unstructured  cells. 

As  mentioned  previously,  the  burning  propellant  surface  boundary  condition  was  used,  along  with  impermeable 
wall  boundaries  to  represent  solid  walls.  The  downstream  exit  boundary  condition  was  an  outflow  with  an  imposed 
back  pressure.  The  computational  domain  was  initialized  with  a  given  pressure  and  temperature  to  start  the  steady 
state  simulation.  The  model  was  executed  in  local  time  stepping  mode  with  a  physical  time  step  of  0.01  seconds. 
This  time  step  was  dynamically  reduced  by  the  Loci-CHEM  solver  in  the  event  that  a  local  Courant  -Friedrichs- 
Lewy  (CFL)  number  exceeded  5,000  or  a  change  in  temperature,  pressure  or  density  of  more  than  5  percent 
occurred  in  one  time  step.  Using  these  settings,  the  simulation  required  about  2,000  time  steps  to  converge  to  a 
steady  state  solution.  Convergence  was  evaluated  by  observing  the  drop  in  residuals  and  by  observing  volume 
probes  of  temperature  and  pressure  at  several  locations  throughout  the  flowfield. 


B.  Baseline  Results 


Figure  2  shows  a  comparison  between  the 
Loci-CHEM  solution,  the  baseline  1-D  solution 
(Jacobs’  1-D  analytical  model1),  and  the 
experimental  measurement1  along  with  10 
percent  error  bars  for  debond  angle  and 
pressure  chamber  of  4.75°  and  400  psi 
respectively.  The  results  show  excellent 
agreement  among  all  three  results.  The  1-D 
model  and  3-D  CFD  results  are  virtually 
identical  and  both  profiles  pass  neatly  through 
the  measurement  error  bars  (downstream 
measurement  results  were  not  available  for 
comparison).  There  is  a  small  deviation 
between  the  fluid  dynamic  models  near  the 
crack  tip,  as  the  1-D  model  predicts  a  debond 
tip  pressure  of  about  500  psi,  whereas  the  3-D  Figure^.  Comparison  between  analyses  and  experiment 
CFD  solution  predicts  about  485  psi  at  the  tip  at  4-75  and  400  psi. 

(-2.5%  difference,  which  is  likely  within  the 

uncertainty  of  both  the  CFD  solution  and  1-D  model).  Experimental  data  were  not  collected  at  the  debond  tip,  so  it 
is  not  apparent  which  solution  is  correct.  Similar  pressure  profile  agreement  is  also  observed  for  the  same  debond 

angle  with  chamber  pressures  of  300  and  170  psi. 


Distance  along  debond  -  inch 


Distance  along  debond  -  inch 

Figure  3.  Comparison  between  analyses  and  experiment 
at  3.17°  and  330  psi. 

a  difference  of  -18%  in  pressure  at  the  debond  tip. 
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Figure  3  shows  the  results  comparison  for  a 
smaller  debond  angle  and  lower  chamber 
pressure.  For  this  analysis,  the  debond  angle  was 
reduced  to  3.17°  and  the  chamber  pressure  was  set 
to  about  332  psi.  There  is  slightly  a  higher 
difference  of  the  results  of  the  pressure  calculation 
at  the  debond  tip  of  -6%,  which  might  be  caused 
by  a  higher  Mach  number  flow  inside  the  debond 
due  to  a  smaller  volume.  However,  an  overall 
good  agreement  was  obtained  for  this  analysis. 

One  particular  experiment  was  purposely 
designed  to  induce  choked  flow  within  the 
debonded  volume  and,  as  a  result,  the  fluid 
dynamics  changed  significantly.  Figure  4  shows 
the  comparisons.  In  this  instance,  the  CFD 
solution  appears  to  be  better,  compared  to  the 
experimental  data,  than  the  1  -D  solution.  There  is 
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There  were  several  important  conclusions  to  this  baseline  set  of  analyses.  First,  the  CFD  represented  the  general 
physics  of  the  particular  problem  faithfully:  the  difference  in  bulk  flowfield  behavior  between  low  subsonic  flow 
(high  wedge  angle,  low  chamber  pressure)  and  sonic  flow  (0°  wedge  angle,  high  chamber  pressure)  is  clearly 
captured  by  the  CFD  model  with  no  special  additional  input  required.  Second,  despite  some  small  quantitative 
differences  between  experimental  measurements 
and  CFD  results,  the  correspondence  between  the 
two  is  qualitatively  excellent.  To  better 
understand  those  quantitative  differences,  a  subset 
of  analyses  was  completed  that  highlighted  the 
sensitivity  of  the  CFD  results  to  several 
assumptions. 

As  discussed  previously,  the  flowfield  inside 
the  simulated  flaw  is  sensitive  to  the  actual  size  of 
the  crack  tip  gap  and  to  the  overall  empty  volume, 
both  of  which  constantly  changed  throughout  a 
test  due  to  the  propellant  regression. 

Unfortunately,  the  researchers  did  not  record  in 
Ref.  1  exactly  when  any  particular  pressure 
measurement  was  made  in  any  given  test.  They  Figure  4.  Comparison  between  analyses  and  experiment 

merely  stated  that  each  dataset  came  from  within  a(  q°  an(j  245  psi. 

the  first  0. 1  to  0.5  seconds  of  any  given  test.  This 

made  it  impossible  to  determine  whether  the  idealized  numerical  domain,  shown  in  the  figure  below,  truly 
represented  the  flaw  tip  gap  and  overall  empty  volume  for  any  of  the  experimental  results.  In  every  case,  the 
numerical  domain  used  pre-test  dimensions  of  the  hardware  and  propellant  slab. 

However,  it  was  also  important  to  understand  better  how  several  assumptions  affected  the  CFD  results,  since 
future  analyses  on  real  motor  flaws  will  have  analysts  grapple  with  much  of  the  same  geometric  uncertainty.  In  fact, 
they  will  likely  face  uncertainties  in  flaw  geometry  to  a  far  greater  extent,  since  NDE  equipment  generally  do  not 
have  high  resolution. 

Therefore,  several  sensitivity  analyses  were  completed.  First,  a  grid  sensitivity  study  was  performed  to 
understand  whether  grid  independence  had  been  achieved.  This  was  followed  with  an  analysis  to  determine  how 
sensitive  the  results  were  to  assumptions  of  gas  properties.  Finally,  the  effect  of  turbulence  level  assumptions  was 
analyzed.  All  of  these  results  are  shown  in  the  following  sections. 


C.  Grid  Sensitivity  Study 

Figure  5  shows  a  comparison  between  two 
different  grid  densities;  the  As  term  is  the  factor  that 
controls  grid  spacing.  In  this  instance,  the  spacing 
factors  of  0.01  and  0.005  are  used  (i.e.,  the  second  grid 
used  spacing  half  as  much  as  the  first),  which  roughly 
translates  into  a  mesh  density  increase  by  a  factor  of 
eight,  via  a  doubling  of  the  number  of  grid  points 
along  all  three  primary  dimensions.  The  results 
between  the  two  analyses  are  nearly  identical,  showing 
that  a  grid  independent  solution  was  obtained. 

D.  Combustion  Gas  Properties  Sensitivity  Study 

Even  with  very  complex  solid  propellant 
simulations  addressing  other  difficult  problems  facing 
solid  rocket  motors,  complex  chemical  kinetic  models 
are  not  often  necessary7.  It  was  assumed  that  such  will 
be  the  case  for  the  flaw  propagation  problem,  since  gross  flowfield  properties  of  pressure  and  Mach  number  are  only 
moderately  affected  by  chemical  non-equilibrium  (local  temperature  is,  of  course,  the  primary  affected  property). 
Thus  the  baseline  results  reported  here  reflect  the  assumption  of  chemical  equilibrium.  Properties  were  generated 
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using  the  NASA  CEA2  code,  the  current  version  of  the  venerable  code  developed  by  Gordon  and  McBride  starting 
in  the  early  1960’s.  CEA2  was  used  to  calculate  specific  heat  (Cp),  ratio  of  specific  heats  (y),  adiabatic  flame 
temperature  and  density.  Recall  also  that  there  was  considerable  uncertainty  in  the  actual  propellant  tested  in  Ref.  1. 


To  better  understand  how  sensitive  the  baseline 
results  were,  an  additional  analysis  was  completed 
that  used  the  (frozen)  combustion  gas  properties 
reported  in  Ref  1 .  Figure  6  shows  that  it  does  make  a 
small  difference  to  the  CFD  results,  especially  near 
the  flaw  tip  region,  but  all  results  fell  within  the 
estimated  uncertainty  of  the  measurements.  For  the 
baseline  analysis,  it  appears  that  equilibrium 
assumptions  are  satisfactory. 

E.  Turbulence  Level  Sensitivity  Study 

After  an  extensive  review,  the  literature  does  not 
show  a  discernible  trend  regarding  turbulence 
assumptions  for  burning  flaw  flowfield  simulation. 
In  fact,  the  review  shows  a  roughly  even  division 
between  simulations  that  assumed  laminar  flow  and 
ones  that  assumed  turbulent  flow.  Given  that  solid 
rocket  propulsion  typically  exhibits  highly  turbulent 
flow  within  the  grain  perforation,  it  seems  natural  to 
assume  that  the  flow  within  a  flaw  would  likewise  be 
turbulent.  However,  it  also  seems  reasonable  to 
assume  that  the  flow  will  transition  from  laminar  to 
turbulent,  similar  to  any  non-combusting  channel  or 
pipe  flow. 

To  understand  the  sensitivity  of  the  CFD  results 
to  turbulence  assumptions,  an  analysis  was  completed 
using  a  fully  laminar  assumption.  The  results,  shown 
in  Figure  7,  show  a  higher  pressure  near  the  crack  tip, 
due  to  lower  irreversible  losses;  as  with  the  other 
sensitivity  analyses,  the  difference  seen  is  well  within 
the  experimental  uncertainty.  The  best  understanding 
at  the  moment  is  that  that  the  flow  inside  a  burning 
flaw  will  probably  transition  from  laminar  to  turbulent. 


Figure  6.  Combustion  gas  properties  sensitivity 
study  results. 


Distance  along  debond  -  inch 

Figure  7.  Turbulence  intensity  sensitivity  study 
results. 


V.  Conclusions 

The  goal  of  the  current  study  was  to  demonstrate  feasibility  of  using  a  modern  CFD  tool  to  predict  the  flowfield 
inside  a  solid  propellant  grain’s  burning  flaw.  This  is  the  first  step  to  realizing  the  use  of  a  multiphysics  approach  to 
simulate  flawed  solid  rocket.  The  baseline  analyses,  along  with  the  understanding  of  the  drivers  behind  the  minor 
discrepancies  shown  between  the  CFD  results  and  the  experimental  results,  clearly  demonstrate  feasibility  of  such 
approach.  Sensitivity  analyses  showed  that  assumptions  regarding  grid  resolution,  chemical  non-equilibrium,  and 
turbulence  level  are  relatively  modest  for  such  a  simple  geometry.  Actual  flaws  are  not  expected  to  show  low 
sensitivity.  Realistically,  specification  of  turbulence  levels  likely  provide  the  largest  source  of  sensitivity,  so  work 
remains  in  understanding  the  stability  of  the  flow  inside  a  burning  flaw.  Prediction  of  crack  tip  pressure  is  affected 
by  the  level  of  turbulence,  so  assumptions  must  be  accurate  in  order  to  increase  the  accuracy  of  future  flaw 
propagation  predictions. 
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